
refract_time=50;
l_win=1000;
r_win=10000;
min_dur=10;

do_layer=2;

cell_cnt=0;


snip_cnt=0;
snipS_cnt=0;
snipE_cnt=0;
vel_snip_cnt=0;
vel_snip_c_cnt=0;
mo_snip_cnt=0;
vs_snip_cnt=0;
vs_stop_snip_cnt=0;
%act_ind=[];
cell_id=[];
active_cells=[];
xcorr_delay=[];



ps_snips=zeros(l_win+r_win+1,1598);
ps2_snips=zeros(l_win+r_win+1,1598);
co_start_snips=zeros(l_win+r_win+1,1598);
co_stop_snips=zeros(l_win+r_win+1,1598);
mo_start_snips=zeros(l_win+r_win+1,1598);
co_start_fb_snips{1}=zeros(l_win+r_win+1,1598);
co_start_fb_snips{2}=zeros(l_win+r_win+1,1598);
co_start_fb_snips{3}=zeros(l_win+r_win+1,1598);
% 
% 
%     mo_snips=zeros(l_win+r_win+1,150000);
%     vs_snips=zeros(l_win+r_win+1,250000);
%     vs_stop_snips=zeros(l_win+r_win+1,100000);
%     saccade_snips=zeros(l_win+r_win+1,300000);


ps_snip_time=[];
vel_pb_session=[];


pcf=[];



for ind=1:length(AV_data)
    disp(num2str(ind))
    if AV_data(ind).layer==do_layer & AV_data(ind).ACode==1
        
        type=(AV_data(ind).velP_ind&AV_data(ind).velM_ind)+2*(AV_data(ind).velP_ind&~AV_data(ind).velM_ind)+3*(~AV_data(ind).velP_ind&AV_data(ind).velM_ind);
        
        conc_offsets=intersect(find(diff(AV_data(ind).velP_ind)==-1),find(diff(AV_data(ind).velM_ind)==-1));
        co_stop=10*conc_offsets(logical([~diff(conc_offsets)'<refract_time 1]));
        
        conc_onsets=intersect(find(diff(AV_data(ind).velP_ind)==1),find(diff(AV_data(ind).velM_ind)==1));
        conc_onsets=conc_onsets(conc_onsets>refract_time);
        use_onsets=logical([]);
        for knd=1:length(conc_onsets)
            use_onsets(knd)= sum(AV_data(ind).velP_ind(conc_onsets(knd)-refract_time:conc_onsets(knd)))==0 & sum(AV_data(ind).velM_ind(conc_onsets(knd)-refract_time:conc_onsets(knd)))==0;
        end
        co_start=10*conc_onsets(use_onsets);
        
        ps_raw=find_ps_during_running(AV_data(ind).velP_ind,AV_data(ind).velM_ind);
        ps_start=parse_cond_starts(ps_raw,min_dur,refract_time,1);
        
%         ps_start=find(diff((AV_data(ind).PS)>0.75)==1);
%         ps_start=ps_start(AV_data(ind).velP_ind(round(ps_start/10)));
        
        ps2_start=parse_cond_starts(logical(AV_data(ind).velP_ind&~AV_data(ind).velM_ind - ps_raw),min_dur,refract_time,1);
        
        mo_start=parse_cond_starts_vs(AV_data(ind).velP_ind,min_dur,refract_time,1);
        
        vs_start=parse_cond_starts_vs(~AV_data(ind).velP_ind&AV_data(ind).velM_ind,min_dur,refract_time,1);
        
        vs_raw=find_pb_during_sitting(AV_data(ind).velP_ind,AV_data(ind).velM_ind);
        vs_stop=parse_cond_starts_vs(vs_raw,min_dur,refract_time,2);
        
        
        pb_onsets=[];
        if length(AV_data(ind).shutter_open)==2
            pb_onsets(1)=AV_data(ind).shutter_open(1);
            pb_onsets(2:4)=find(diff(AV_data(ind).PS(AV_data(ind).shutter_open(2):end)>0.1)>0)+AV_data(ind).shutter_open(2);
            pb_offsets=find(diff(AV_data(ind).PS(AV_data(ind).shutter_open(2)+1000:end)>0.1)<0)+AV_data(ind).shutter_open(2);
            pb_dur=pb_offsets(1)-pb_onsets(2);
        elseif length(AV_data(ind).shutter_open)==3
            pb_onsets(1)=AV_data(ind).shutter_open(1);
            pb_onsets(2:4)=find(diff(AV_data(ind).PS(AV_data(ind).shutter_open(3):end)>0.1)>0)+AV_data(ind).shutter_open(3);
            pb_offsets=find(diff(AV_data(ind).PS(AV_data(ind).shutter_open(3)+1000:end)>0.1)<0)+AV_data(ind).shutter_open(3);
            pb_dur=pb_offsets(1)-pb_onsets(2);
        elseif length(AV_data(ind).shutter_open)>=4
            pb_onsets=AV_data(ind).shutter_open+1000;
            pb_dur=pb_onsets(2)-pb_onsets(1)-2000;
        else
            pb_onsets=[];
        end
        
        if length(AV_data(ind).shutter_open)==2
            tmp_PS=round(3*AV_data(ind).PS);
            tmp_PS(1:AV_data(ind).shutter_open(2))=1;
            fb_onsets=find(diff(tmp_PS==0)==1);
            fb_offsets=find(diff(tmp_PS==0)==-1);
        else
            fb_onsets=[];
            fb_offsets=[];
        end

        ps_start=ps_start(ps_start<AV_data(ind).shutter_open(2));
        ps2_start=ps2_start(ps2_start>AV_data(ind).shutter_open(2));
        
        use_ps=zeros(length(ps2_start),1);
        for knd=1:length(ps2_start)
            if sum(AV_data(ind).velP_ind(round(ps2_start(knd)/10)-500:round(ps2_start(knd)/10)-1))==0&sum(AV_data(ind).velM_ind(round(ps2_start(knd)/10)-500:round(ps2_start(knd)/10)-1))==0
                use_ps(knd)=1;
            end
        end
        ps2_start=ps2_start(logical(use_ps));
        
        
        
        for knd=1:size(AV_data(ind).raw_act,2)
            if max(AV_data(ind).raw_act(:,knd))<40
                cell_cnt=cell_cnt+1;
                std_noise=sqrt(sum((AV_data(ind).raw_act(AV_data(ind).raw_act(:,knd)<1,knd)-1).^2/sum(AV_data(ind).raw_act(:,knd)<1)));
                
                [ROI_act_ms]=map_ROI_act_to_aux_data(AV_data(ind).raw_act(:,knd),AV_data(ind).frame_times,length(AV_data(ind).PS));
                ps_cnt_tmp=0;
                for jnd=1:length(ps_start)
                    try
                        if ~sum((ROI_act_ms(ps_start(jnd)-l_win:ps_start(jnd)+r_win))==0)
                            ps_cnt_tmp=ps_cnt_tmp+1;
                            ps_snips(:,cell_cnt)=ps_snips(:,cell_cnt)+ROI_act_ms(ps_start(jnd)-l_win:ps_start(jnd)+r_win);
                            
%                             ps_snip_time(snip_cnt,1)=ps_start(jnd);
%                             ps_snip_time(snip_cnt,2)=jnd;
%                             pb_session=find(sign(ps_start(jnd)-pb_onsets)>0,1,'last');
%                             ps_snip_time(snip_cnt,3)=ps_start(jnd)-pb_onsets(pb_session);
%                             ps_snip_time(snip_cnt,4)=pb_session;
%                             ps_snip_time(snip_cnt,5)=mean(AV_data(ind).velP(ps_start(jnd)-500:ps_start(jnd)));
%                             ps_snip_time(snip_cnt,6)=interim_runnig(pb_session);
%                             ps_snip_time(snip_cnt,7)=ps_start(jnd)-mo_start(find(ps_start(jnd)-mo_start>0,1,'last'));
%                             ps_snip_time(snip_cnt,8)=mean(AV_data(ind).velP(ps_start(jnd):ps_start(jnd)+3000));
%                             ps_snip_time(snip_cnt,9)=ind;
%                             ps_snip_time(snip_cnt,10)=mean(AV_data(ind).velP(ps_start(jnd)-400:ps_start(jnd)-200));
%                             ps_snip_time(snip_cnt,11)=mean(AV_data(ind).velM(ps_start(jnd)-400:ps_start(jnd)-200));
%                             ps_snip_time(snip_cnt,12)=mean(AV_data(ind).velP(ps_start(jnd):ps_start(jnd)+200));
%                             ps_snip_time(snip_cnt,13)=mean(AV_data(ind).velP(ps_start(jnd)+501:ps_start(jnd)+1500));
%                             ps_snip_time(snip_cnt,13)=cell_cnt;
%                             ps_snip_time(snip_cnt,14)=knd;
                            
                        end
                    end
                end
                ps_snips(:,cell_cnt)=ps_snips(:,cell_cnt)/ps_cnt_tmp;
                
                
                ps2_cnt_tmp=0;
                for jnd=1:length(ps2_start)
                    try
                        if ~sum((ROI_act_ms(ps2_start(jnd)-l_win:ps2_start(jnd)+r_win))==0)
                            ps2_cnt_tmp=ps2_cnt_tmp+1;
                            ps2_snips(:,cell_cnt)=ps2_snips(:,cell_cnt)+ROI_act_ms(ps2_start(jnd)-l_win:ps2_start(jnd)+r_win);
                        end
                    end
                end
                ps2_snips(:,cell_cnt)=ps2_snips(:,cell_cnt)/ps2_cnt_tmp;
                
                nbr_ps2(cell_cnt)=length(ps2_start);
                
                co_start_cnt_tmp=0;
                co_start_fb_cnt_tmp=[0 0 0];
                for jnd=1:length(co_start)
                    try
                        if ~sum((ROI_act_ms(co_start(jnd)-l_win:co_start(jnd)+r_win))==0)
                            co_start_cnt_tmp=co_start_cnt_tmp+1;
                            co_start_snips(:,cell_cnt)=co_start_snips(:,cell_cnt)+ROI_act_ms(co_start(jnd)-l_win:co_start(jnd)+r_win);
                                                        
                            curr_fb_session=find(co_start(jnd)>fb_onsets&co_start(jnd)<fb_offsets);
                            if isempty(curr_fb_session)
                                curr_fb_session=0;
                            end

                            if curr_fb_session>0
                                co_start_fb_cnt_tmp(curr_fb_session)=co_start_fb_cnt_tmp(curr_fb_session)+1;
                                co_start_fb_snips{curr_fb_session}(:,cell_cnt)=co_start_fb_snips{curr_fb_session}(:,cell_cnt)+ROI_act_ms(co_start(jnd)-l_win:co_start(jnd)+r_win);
                            end
                        end
                    end
                end
                co_start_snips(:,cell_cnt)=co_start_snips(:,cell_cnt)/co_start_cnt_tmp;
                for jnd=1:3
                    co_start_fb_snips{jnd}(:,cell_cnt)=co_start_fb_snips{jnd}(:,cell_cnt)/co_start_fb_cnt_tmp(jnd);
                end
                
                co_stop_cnt_tmp=0;
                for jnd=1:length(co_stop)
                    try
                        if ~sum((ROI_act_ms(co_stop(jnd)-l_win:co_stop(jnd)+r_win))==0)
                            co_stop_cnt_tmp=co_stop_cnt_tmp+1;
                            co_stop_snips(:,cell_cnt)=co_stop_snips(:,cell_cnt)+ROI_act_ms(co_stop(jnd)-l_win:co_stop(jnd)+r_win);
                        end
                    end
                end
                co_stop_snips(:,cell_cnt)=co_stop_snips(:,cell_cnt)/co_stop_cnt_tmp;
                
                mo_start_cnt_tmp=0;
                
                
                for jnd=1:length(mo_start)
                    try
                        if ~sum((ROI_act_ms(mo_start(jnd)-l_win:mo_start(jnd)+r_win))==0)
                            mo_start_cnt_tmp=mo_start_cnt_tmp+1;
                            mo_start_snips(:,cell_cnt)=mo_start_snips(:,cell_cnt)+ROI_act_ms(mo_start(jnd)-l_win:mo_start(jnd)+r_win);

                        end
                    end
                end
                
                mo_start_snips(:,cell_cnt)=mo_start_snips(:,cell_cnt)/mo_start_cnt_tmp;

             
%                 act=AV_data(ind).raw_act(:,knd)-1;
%                 act=act.*(act>(3.72*std_noise))+1;
%                 
                for jnd=1:4
                    cell_id(1,cell_cnt)=ind;
                    cell_id(2,cell_cnt)=knd;
                    cell_id(3,cell_cnt)=AV_data(ind).ExpID;
                    cell_id(4,cell_cnt)=AV_data(ind).same_site_as;
                end

            else
                disp(['activity in cell ' num2str(AV_data(ind).ExpID) ' - ' num2str(knd) ' too large'])
            end
        end
        
        
    end
end


bad_cells=max(abs(ps_snips))>6;

ps2_snips(:,bad_cells)=nan;
ps_snips(:,bad_cells)=nan;
mo_start_snips(:,bad_cells)=nan;
co_stop_snips(:,bad_cells)=nan;
co_start_snips(:,bad_cells)=nan;

for ind=1:3
    co_start_fb_snips{ind}(:,min(co_start_fb_snips{ind})<0.1)=nan;
end

% ps_fb_mean=nanmean(ps_snips(:,nns)');
% ps2_fb_mean=nanmean(ps2_snips(:,nns([1:13 15:16 18:31]))');
% 
% 
% ps_fb_mean=nanmean(ps_snips');
% ps2_fb_mean=nanmean(ps2_snips');
% 
% figure;
% plot(ps_fb_mean-mean(ps_fb_mean(501:1000)))
% hold on
% plot(ps2_fb_mean-mean(ps2_fb_mean(501:1000)))
% 
% xlim([500 6000])
% box off
csfs=[];
for ind=1:3
    csfs(:,ind)=nanmean(co_start_fb_snips{ind}');
end

figure;
clrind='rkb';
hold on
for ind=1:3
    plot(csfs(:,ind)-mean(csfs(501:1000,ind)),'color',clrind(ind));
end
xlim([500 4000])

fbact1=mean(co_start_fb_snips{1}(1500:2500,:))-mean(co_start_fb_snips{1}(501:1000,:));
fbact3=mean(co_start_fb_snips{3}(1500:2500,:))-mean(co_start_fb_snips{3}(501:1000,:));

fbact1=fbact1(~isnan(fbact3));
fbact3=fbact3(~isnan(fbact3));


figure;
plot(fbact1,fbact3,'.')

xlim([-0.2 1.2])
ylim([-0.2 1.2])
hold on
plot([-0.2 1.2],polyval(polyfit(fbact1,fbact3,1),[-0.2 1.2]))




ps_mean=nanmean(ps_snips');
mo_start_mean=nanmean(mo_start_snips');
co_stop_mean=nanmean(co_stop_snips');
co_start_mean=nanmean(co_start_snips');

figure;
plot(mo_start_mean-mean(mo_start_mean(501:1000)));
hold on
plot(co_stop_mean-mean(co_stop_mean(501:1000)));
plot(ps_mean-mean(ps_mean(501:1000)),'r');
plot(co_start_mean-mean(co_start_mean(501:1000)),'k');


figure;
plot(nanmean(ps2_snips(:,nns)')-mean(nanmean(ps2_snips(501:1000,nns)')))
hold on
plot(nanmean(ps_snips(:,nns)')-mean(nanmean(ps_snips(501:1000,nns)')),'r')

xlim([1 10000])






